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ABSTRACT 

We investigate the phase space structure of the relativistic Sitnikov problem in the first 
post- Newtonian approximation. The phase space portraits show a strong dependence 
on the gravitational radius which describes the strength of the relativistic pericentre 
advance. Bifurcations appearing at increasing the gravitational radius are presented. 
Transient chaotic behavior related to escapes from the primaries are also studied. 
Finally, the numerically determined chaotic saddle is investigated in the context of 
hyperbolic and non-hyperbolic dynamics as a function of the gravitational radius. 
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1 INTRODUCTION 



The Sitnikov problem (SP) (|Sitnikovll960l ) is one of the sim- 
plest dynamical systems in celestial mechanics that provides 
all kinds of chaotic behavior. The configuration of the system 
is defined by: two point-like bodies of equal masses (called 
primaries) orbiting around their common centre of mass due 
to their mutual gravitational forces, and a third body of neg- 
ligible mass moving along a line, perpendicular to the orbital 
plane of the primaries, going through their barycentre. For 
the circula r motion of the p rimaries, the problem is inte- 
grable and iMacmillanI (|l91lh gave a closed form analytical 
solution with elliptic integrals . iMoseJ l|l973l l showed the ex- 
istence of chaotic behavior using symbolic dynamics. 

In the last decades the problem was i nvestigated 
in de tails both analytically and numerically. iLiu fc SurJ 
derived a ma pping model to investigate the prob- 
Wodnaij (|l99ll ) introduced a new formulation for 
equation of motion by using the t rue anomaly of 
as independent variable. If fagel fc Lhotkal 
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lem. 
the 

the primaries 



l|2005l ) extended the analytical approximations up to very 
high o rders by using extensive computer algebra. iDvorakI 
1 19931 ) showed by numerical computations that invariant 
curv es exist for sma ll oscil lations around the barycen- 
tre. lAlfaro fc Chiraltl (|l993l ) determined invariant rota- 
tional curves by applying the Birkhoff normal form of 
an are a- preserving mapping. Periodic solutions were stud- 



Belbruno etHI |l994) . 



Kallrath et al. ||199' 



ied by Perdios fc MarkcUoa 
Jalali fc Pourtakdoust 41997 

Corbera fc Llibre (2000|). The complete phase space was 



and 



studie d numerically by iDvorakI (|2007l ) and iKovacs fc Erdi I 
l|2007l ). 

The Sitnikov problem is the perfect manifestation of 
a scattering process in which a particle approaches a dy- 
namical system from infinity, interacts with the system and 
ultimately the particle leaves i t. Escapes to in fi nity i n Sit- 
nikov problem were studied by iKovacs fc Erdil (|2009l ). The 
test particle can escape the system via different exits, thus 
one can determine the basins of escape, since in these sys- 
tems infinity acts as an attractor for an escaping particle. 
A detailed s t udy a bout basins of es c ape c an be found in 
iBleher et all ([198^) and IContopoulosI (|2002l ). 

An interesting question is how the structure of the 
phase spa c e ch anges due to relativistic effects. Since 
iRobertsonI (|l938t ) gave the solution of the relativistic two- 
body problem in the post-Newtonian (PN) approxima- 
tion, many papers have dealt with pericentre advance 
in celestial mechanics, especially in the case of binary 
pulsars where t he masses of celestial bodies are of the 
same mag nitude (|Damour fc Schaffeill98^ : IWagoner fc Will 
Il976l ). It is an established fact that the pericentre pre- 
cession in the post-Newtonian tw o-body problem is the 
same as in Schwarzschild 's metric (|Landau fc Lifsij Il975l : 
iDamour fc Deruell^ 1 19851 ) expressed by the total mass. 
Namely, the pericentre advance in one revolution is 



A<j!> = 67r 



fc(mi -I- m2) 



(1) 
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where k represents the gravitational constant, mi and m2 
are the masses of the bodies, and a and e describe the clas- 
sical semi-major axis and the eccentricity, respectively. In 
Eq. ([T} c denotes the speed of light. 

The aim of the present work is to show the relativistic 
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dynamics of the Sitnikov problem by taking into account 
the leading post-Newtonian "perturbation". The paper is 
organized as follows. In Sec. 2 we describe the model: first, 
the post-Newtonian two-body problem, then the relativistic 
Sitnikov problem. Section 3 contains our numerical results. 
Sec. 3.1 concentrates on the phase space structure of the rel- 
ativistic SP, and 3.2 deals with chaotic scattering. In Section 
4 we draw our conclusions. The derivation of the equations 
of motion can be found in the Appendix. 



2 DESCRIPTION OF THE MODEL 

The Sitnikov problem is a particular case of the restricted 
three-body problem. The third, mass-less body has no effect 
on the primaries' motion, neither in the classical, nor in the 
relativistic case. Therefore, the problem can be split into 
two parts. The solution of the two-body problem is needed 
for the determination of the motion of the test-particle. The 
instantaneous position of primaries provides the time depen- 
dent driving acting on the test-particle. 

2.1 Post-Newtonian two-body problem 

All our results are given in the leading post-Newtonian (PN) 
approximation based on an assumption of weak inter-body 
gravitational field and slow orbital motions. Beyond the clas- 
sical limit, it contains terms of order v'^ /c? , wh ere t> is a typ- 
ical o rbital velocity and c is the speed of light jCalura et al.l 
Il997f ). 

Throughout the paper the length unit will be cho- 
sen as the semi-major axis of the classical two-body prob- 
lem: a = —kinim2/{2Ec), where Ec is the classical energy 
in the centre of mass frame. The time unit is taken as 
T = a^^'^[k{mi + 7712)]"^/^ from Kepler's third law. The en- 
ergy unit will be km\m2/a. The total dimensionless classical 
energy becomes Ec — —0.5. 

The equations of the relative motion in the post- 
Newtonian centre o f mass frame can be deriv ed from the 
Lagrangian given in iDamour fc Deruelld l|l985l ). The equa- 
tion for the relative coordinate r reads in dimensionless form 
up to first order in A 



r = V = 7 + A 



.(l + 3^)_^^ + _^_(rv) 
+ {4 + 2u)^ + {4-2u)^{rV 



(2) 



where u = mim2 / (nii + 1x12)^ is the effective mass and 

X = k{mi + m2) / ac^ (3) 

denotes the dimensionless gravitational radius with a as the 
classical semi-major axis. 

The invariance of the Lagrangian under time translation 
and spatial rotation implies the existence of four first inte- 
grals, the energy and the angular momentum of the binary 
system in the centre of mass frame: 



E =-v 
2 



1 



3u)v 



+ 



in , \ 2 , (rv) 



(4) 



J =r X V 



1 + A( i(l-3//)7;2 + (3 + j^)i 
' 2 r 



This implies that the motion reduces to an effectively one- 
dimensional bounded problem that in turn leads to a pe- 
riodic time dependence for all the variables. One can also 
derive the dimensionless form of the period P between tw o 
consecutive pericentre passages (|Damour fc Schaffeilll988l ). 
The pericentre advance and P are given as 



A 



27r 



1 - 



(-2i;)3/2 



1 + 7(15 



A 



(5) 



The classical eccentricity is obtained from Ec and the clas- 
sical angular momentum Jc as 



ec = 1 - J< 



21 1/2 



(6) 



The relativistic orb ital elements are then 
iDamour fc SchaffeillTgsi ) 



«W=-^-(^-7)^A, 



e(A) 



1 + 2E.f - A (i/ - 6) 



15 
2 



1/2 



(7) 



Since the masses of primaries are equal, we can set the value 
oi V = 1/4. The solution of Eqs. ([2]) provides the time de- 
pendent driving for the relativistic Sitnikov problem. 



2.2 The relativistic Sitnikov problem 

To obtain the equation of the relativistic Sitnikov problem 
(RSP) we use the Lagrangian o f the post-Newtonian 3-body 
system (|Landau fc Lifsidll975h . The final form of the equa- 
tion is (for details see the Appendix.): 



+ A 



16 z 



4 rp '-' ^ 2 p 



-H6- 



3 (vr)i 
2 p3 



3 (vr)2 



16 



(8) 



here z and i are the dimensionless position and the velocity 
of the test particle, respectively, and p = \J z^ -\- r"^ j^. The 
additional terms in the bracket on the right hand site of 
Eq. ((HI) describe the relativistic effects on Sitnikov's motion. 



3 NUMERICAL RESULTS 

As initial conditions to the two-body problem we take the 
initial conditions at the pericentre in the form rperic = 
(1 - ec; 0), V = (^0; ((1 + ec)/(l - ec)f''^^ ■ We shaU fix the 
classical eccentricity to be ec = 0.2 (Jc = 0.9797). The total 
energy is then 



1 A 71 58ec 3eg 
2'^32 (l-ec)2 



-0.5-^4.04A. 



The period P up to first order in A is therefore 
P = 27r(l -f 14A). 



(9) 



(10) 



Our numerical investigations show that the post- 
Newtonian approximation is valid between and some Ac. 
We define the critical Ac as a value where the numerical 
results of the two-body problem differ from the analytical 
results ((5)1 by about 10 per cent. Although the equations of 
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Figure 1. Analytical and numerical solutions of the relativistic 
two-body problem. Upper panel: Period P. solid line — numeri- 
cal solution, dash-dotted — analytical solution given by Eq. mOI I. 
Lower panel: Pericentre advance, solid line - numerical solution, 
dash-dotted — analytical solution given by Eq. l|5)l. 
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Figure 2. Stroboscopic phase portrait of the RSP (A = 0.005) 
taken with the classical period of the two-body problem. Initial 
conditions are taken from the interval 0.05 ^ z ^ 3.0 (Az = 0.05) 
and 2 = 0. 



motion ((2} are valid up to first order in A, they are nonlin- 
ear equations and can provide results which are higher order 
in A or /(?. There is thus a threshold beyond which the 
numerical solution no longer holds. Figure [1] indicates that 
Ae « 0.035. 

Due to the time periodic driving one can investigate the 
structure of the three dimensional phase space via strobo- 
scopic or Poincare maps, like in the non-relativistic case. 



3.1 Phase space structure 

Along with the gravitational radius the orbital period 
changes, as expressed by Eq. (llOf) . Therefore, if data are 
stored corresponding to the Keplerian orbital period, we ob- 
tain a confused phase portrait. Looking at Figure [2J one 
cannot distinguish islands or chaotic bands, we see just 
"fuzzy" curves and sparse points everywhere. In order to get 
a "transparent" phase portrait, the new period, P, of the pri- 
maries' revolution or a redefined Poincare map is needed. In 
this paper figures are plotted corresponding to the Poincare 
map taken at r = Vperic, when the primaries are in the peri- 
centre. 

Figure [3] shows this Poincare map that allows us to in- 
vestigate the phase space as usual. One can indeed see the 
pattern typical of conservative dynamicsQ By comparing the 
relativistic (Fig.|3} and the classical phase portraits (Fig.|4]), 
the structure is different but the Hamiltonian characteristics 
remain. 

The RSP has one new parameter, the gravitational ra- 
dius A. The qualitative features of the phase space depends 
then on A only. One can see from Fig. [3] and |3] that even 
a small A can change the phase portrait dramatically. For 
A = 0.005 the main difference to the classical case is the 
modified surroundings of the 2:1 resonance. The central re- 
gion is similar to the classical case, invariant tori are situated 




Figure 3. Phase portrait of the RSP at A = 0.005 taken at the 
pericentre passage of the primaries. Initial conditions are the same 
as in Fig. |2] 
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Figure 4. Phase space portrait of the classical Sitnikov problem 
at e = Ec = 0.2, taken at times 2it , An , 6n , . . . Initial conditions 
are the same as in Fig. [2] 



^ For simplicity, instead of the canonically conjugated {z,pz) co- 
ordinates we use the traditional (z, z) coordinates, although it 
makes the Poincare map not exactly area preserving. This map 
is smoothly conjugated to the area preserving one. 



4 T. Kovdcs, Gy. Bene, and T. Tel 



around the stable origin. However, the main islands corre- 
sponding to the 2:1 resonance are changed. The irregularly 
scattered points around the islands of stability represent the 
trajectories that can escape the system sooner or later. 

Figure [5] shows three phase portraits for different values 
of the gravitational radius A. One can draw several conclu- 
sions from these plots. First, it is evident that the central 
stable region is shrinking when A becomes larger. It means, if 
the perturbation is larger, the domain of the ordered motion 
is smaller. Second, the size of the island of the 2:1 resonance 
becomes larger along the direction parallel to the z axis and 
smaller perpendicular to it. During this process, the perime- 
ter of the island increases resulting in a longer island chain 
around the last KAM-torus. Consequently, there are longer 
Cantori which the scattered trajectories can stick to during 
the scattering process. We will see this in the next Section. 

We have identified a bifurcation at A ~ 0.03. The qual- 
itative change due to the increasing gravitational radius is 
well-seen in Figure [B] Panel (a) shows the right island of 
the 2:1 resonance with one stable elliptic fixed point on the 
phase portrait sitting in the middle of the regular island. 
However, for larger A (panel (b) and (c)) the stable periodic 
orbit becomes unstable and two new elliptic fixed points ap- 
pear to the left and right. 

In order to show another picture about how the phase 
space changes with A, we plotted the escape times in a con- 
tour plot (Figure [71 . If the mechanical energy of the test- 
particle becomes positive, the trajectory never returns to the 
primaries' plane, i.e. the test-particle leaves the system. In 
Figure [7] different colors represent different escape times in 
the (A, z) plane. The contour plot was made as follows. We 
have chosen 250 initial conditions along the z-axis in [0; 8] 
with initial velocity zero, and 500 values of A from the inter- 
val [0; 0.035]. The escape times i.e. the time needed to reach 
the state where the energy of the particle becomes positive, 
were computed at the grid points of the (A, z) lattice and 
plotted with different colors. Figure [7] allows us to see the 
evolution of the extent of the stable islands (where the life- 
time is maximal) and filamentary structures. One can see 
a pitchfork bifurcations when A is growing at A ~ 0.032. 
Beyond the value A ~ 0.012 the island of 1:1 resonance 
[z ~ 0.75) becomes separated from the main central stable 
region. In other words, trajectories between the resonant 
island and the central invariant curves may escape to the 
infinity (Fig. [7] and Fig. [S^) . The bright colored filaments 
with higher escape times correspond to a fractal set repre- 
senting the stable manifol d of a chaotic saddle existing far 
from the stability islands (|Kovacs fc Erdil|2009l ). Trajecto- 
ries originating from these initial conditions can spend very 
long time around the primaries' plane before leaving the 
system. 

We can say that the parameter A plays a similar role 
in the RSP as the eccentricity in the classical case. Vary- 
ing the gravitational radius beyond a fixed eccentricity, we 
obtain qualitative changes in the phase space as is com- 
mon in Hamiltonian dynamics. A complete picture of the 
p hase space of th e classical Sitnikov problem was published 
m iDvorakI l|2007l ) where escape times show a structure sim- 
ilar to that of Fig. [7] 




z 



Figure 7. Escape times of the trajectories that originate from 
{z, i = 0) on the (A, z) parameter plane. When the test-particle 
has positive energy, we store the integration time as the escape 
time of the orbit. This contour map was calculated over 100 pe- 
riods of primaries (i.e. roughly 628 time units, see the color bar 
on the right hand side). The central stable region and the right 
island of the 2:1 resonance (around z 2.0) are plotted in light 
gray (orange on-line). These trajectories belong to periodic or 
quasi-periodic orbits and never escape. However, there are other 
regions far from the stable islands where the escape times are 
higher than many orbital periods of primaries. These parts of the 
phase plane contain the stable manifold of the chaotic saddle. 

3.2 Chaotic scattering 

In the previous Section we have seen that there are initial 
conditions which correspond to long life times. Transient 
chaos a ppears as a sca t tering process i n conservative dy- 
namics (lEckhardt|[l987l : ljung fc Schoi3 Il987l : iBleher et al.l 
Il989h . In our example the test particle comes close to the 
primaries' plane and makes several oscillations before escap- 
ing. One can ask where the long-lived trajectories are in the 
phase space. In order to answer this question, a large num- 
ber of points is distributed uniformly in the phase space 
and their evolution in time is followed. We are interested in 
non-escaping trajectories in a preselected region. Before the 
trajectories leave the system they draw out a well-defined 
fractal set in surfaces of sections (|Ot^ll993D . The invari- 
ant object in the phase space responsible for the transient 
chaotic behavior is the chaotic saddle. They characterize 
the dynamics in a way chaotic bands characterize perma- 
nent chaos. It was also shown that these invariant saddles 
have two different parts. One of them, the hyperbolic part, 
is responsible for short lifetimes, the other one, the non- 
hyperbolic part, is situated close to the border of the KAM- 
tori and is associated with the sticky orbits l|Tel fc Gruij 
l2006l : lAhmann fc^l2008l ) . 

3.2.1 Short time escape 

In truly hyperbolic systems the nu mber of non-escaping 
traje ctories decreases exponentially l|Kantz fc GrassbergeJ 
Il985l ). However, in the Sitnikov and also in the relativis- 
tic Sitnikov problem we have regular islands in the phase 
portraits, i.e the phase space is mixed. There are both hy- 
perbolic and non- hyperbolic parts present. In this case the 
decrease of the survivors follows the exponential rule only 
for shorter times, as shown in Figure [8] 
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Figure 5. Phase portraits for different gravitational radii (a):A = 0.015, (b):A = 0.025, {c):A = 0.035. The central stable region becomes 
smaller from left to right. Moreover, the different shape of the island of the 2:1 resonance shows the sensitivity to A. Panel (c) exhibits 
a phase space section after a bifurcation at A = 0.032. 
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Figure 6. Bifurcation due to changing the gravitational radius about A Ki 0.03 Panel (a)-(c): A = 0.028, A = 0.03, A = 0.032. The 
originally stable periodic orbit (2:1 resonance) becomes unstable and two new stable fixed points appear in the phase portrait. 




Table 1. Escape rates and average lifetimes. The greater A the 
shorter the average lifetime of chaos. The function k(A) is close 
to be linear. 



A 


escape rate,/? 


escape time, 1/k 


errors 


0.0 


0.007 


142.8 


±1.4 X 10-5 


0.005 


0.008 


125.0 


±1.5 X 10-5 


0.015 


0.010 


92.6 


±2 X 10-5 


0.025 


0.012 


80.6 


±3.8 X 10-5 


0.035 


0.013 


75.8 


±4.7 X 10-5 



Figure 8. Number N{t) of non-escaped trajectories from a pre- 
selected box (—10 2: ^5 10) in the phase space. Different marks 
represent different values of A : triangles - 0; diamonds - 0.005; 
crosses - 0.015; asterisks - 0.025; squares - 0.035. The escape rates 
can be obtained from the slope of the fitted lines on log-lin plot, 
see the inset on the right. The corresponding escape rates are 
shown in Table [Tl Initial conditions: 6 ^ z ^ 6.8, |2| ^5 0.1. 



If we choose the initial conditions far from regular is- 
lands (e.g. 6 ^ z ^ 6.8 and —0.1 ^ z ^ 0.1), the dynamics 
can be considered to be hyperbolic. Figure [8] shows the num- 
ber of non-escaped trajectories for different gravitational 
radii. One can see that the first segments of the curves follow 
different straight lines in the log-lin plot, N{t) ~ e-"'. The 
escape rate, k, whose inverse tells us the average lifetime of 
chaos, can be obtained from the slope of these lines. The 
inset shows that the slopes are not equal. For various A we 
get different escape rates, see Table [1] 



The chaotic saddle responsible for the finite time chaotic 
motion has a double- fractal structure. One can consider this 
object as the union of all the hyperbolic unstable periodic 
orbits and the i ntersections of their stable and unstable man- 
ifolds (|Tel fc G ruiz 2006). In other words, the scattered test- 
particle jumps randomly on the saddle before leaving it. One 
can see the numerically determined saddles for various grav- 
itational radii in Figure [5] If the number A'o of initial condi- 
tions is large enough, the implemented method (|Tel fc Gruij 
2006) allows us to visualize the saddle itself. We suppose 
that the initial point of a trajectory lies close to the stable 
manifold of the chaotic saddle and we follow the evolution 
of this point forward in time. After some iteration it must 
be in the vicinity of the unstable manifold of the saddle. A 
suitable integration time to can be determined from the av- 
erage lifetime, 1/k, of chaos. We have chosen to ~ (2 — 4)/«:. 
Consequently, the mid-point taken at t « of the tra- 
jectory should be close to the saddle. In order to generate 
the chaotic saddles of Fig[9l we stored the mid-points of the 
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Figure 10. The log-log plot of the number N{t) of non-escaped 
trajectories. Power-law decay holds in the non-hyperbolic part of 
the saddle. Trajectories that can come close to the outer KAM- 
tori may spend very long time around them. The exponent in the 
power-law decay describes the escape rate of these sticky trajec- 
tories. The value of cr is not as sensitive as k but our results show 
that a small fluctuation can be detected related to the size of 
Cantori surrounding the quasi-periodic regions. The parameters 
and notation are the same as in Fig. [Sl 



trajectories that do not escape a preselected box (l^l ^ 10 
\z\ ^ 2) in time to. 

The size of the chaotic saddle decreases in Fig. [9] when 
A is growing and the scenario in Figure[9]is similar to that in 
Fig. [5] The double Cantor structure is dominant but empty 
holes appear at the sites of regular islands. This is the conse- 
quence of the quasi-periodic motion on tori which is perma- 
nent and certainly not chaotic, but can be arbitrarily close 
to the chaotic saddle. 



3.2.2 Long time escape - stickiness 

For longer times the number of non-escaping trajectories 
does not follow the exponential decay. Instead, one observes 
a power-law decay which is slower than the exponential one, 
N{t) ~ t~'^ . It is an established fact that trajectories which 
come close to the outer border of the stability islands may 
stick to them through the debris of the previously destroyed 
KAM curves. A geometrical consequence of the stickiness 
effect in the phase space is the denser saddle structure. In 
other words, one can observe the remnants of previously de- 
stroyed KAM-tori, the so called Cantori, around the stability 
islands, and this part can be identified as the non-hyperbolic 
part of the chaotic saddle. 

Our numerical investigations support the theoretical re- 
sults, namely, exponent a does not depend on the param- 
eters of the dynamical systems. From fitting straight lines 
to the points between 500 and 1500 in Figure IIOI we find 
(J « 3.6 irrespective of A. 



3.2.3 Basins of escape 

In general, when the test-particle escapes the system, i.e. 
never comes back, it approaches infinity. Therefore, although 
in conservative systems there are no attractors, one can con- 
sider infinity as an attractor of the system. In other words. 



beyond the escape energy, infinity behaves as an attractor 
for those trajectories which leave the system. 

Initializing many initial conditions in a fine rectangu- 
lar grid one can identify the basins of escape. Figure 1111 
In the relativistic Sitnikov problem the test-particle may 
leave the system upward or downward from the primaries' 
plane depending on the initial conditions of the trajectory; 
two basins of escape can be identify, the basins of -|-oo and 
— oo, respectively. Figure [TTk shows these basins. The white 
region contains initial conditions that provides the orbits 
escaping the system upward, points marked with dark gray 
color correspond to the orbits leaving the system downward. 
The third part of the phase portrait (hght gray) belongs to 
regular islands the trajectories never escape from. Consider- 
ing only the basins of escape one can see the very complex 
structure of the boundary. The fine-scale structure of the 
boundary shows that it is not a simple curve, rather a frac- 
tal set. Fractality is in general, a result of chaotic motion 
iMcDonald et al.lll985l : [Ottll 19931 ). 

Let us consider the analogy with dissipative systems 
where the stable manifold of hyperbolic unstable fixed points 
provide the boundary of the basin of attraction. Figure [TTb 
shows the stable manifold of the chaotic saddle, i.e. the ini- 
tial conditions of those trajectories which remain for very 
long time in the system. The correspondence is evident be- 
tween panel (a) and (b). We point out that the fractal basin 
boundary of escape in open Hamiltonian systems, which is a 
result of chaotic scattering, is identical with the stable man- 
ifold of t he chaotic saddle r e sponsible for trans ient chaos. 
(See also lTel fc G"rui3 l|2006l ): lErnst et all (|200S ')') 



4 CONCLUSIONS 

In this work we have investigated the relativistic Sitnikov 
problem numerically. The motivation was to show how the 
structure of the phase plane changes under the influence 
of general relativity. The model contains the flrst post- 
Newtonian relativistic corrections, derived from the leading 
order relativistic Lagrangian. Besides the eccentricity of pri- 
maries the gravitational radius is the new parameter of the 
system that is related to the pericentre shift of the two large 
bodies. 

The calculations show that the problem remains a 
driven system but the new driving period corresponds to 
the relativistic orbital period of the binaries. Therefore a 
new stroboscopic or Poincare map was required to visual- 
ize correctly the phase space structure. We found that the 
Poincare section exhibits the well-known picture of Hamil- 
tonian chaos. Moreover, changing the parameter A, the 
phase portrait's structure shows qualitative changes. We 
have pointed out the shrinking of the central regular domain 
when the gravitational radius becomes larger and simultane- 
ously several bifurcations occurred to the 2:1 resonance. In 
order to investigate the chaotic scattering through escapes, 
we have integrated a large number of initial conditions. Two 
different type of escapes, short and long time escape, were 
distinguished to identify the dual structure of the chaotic 
saddle. Short time escapes belong to the hyperbolic part of 
the saddle, the number of non-escaping trajectories from a 
preselected region decreases exponentially. During the scat- 
tering process such trajectories draw out the double fractal 
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Figure 9. Chaotic saddles for different gravitational radii, (a) A = 0.015; (b) A = 0.025; (c) A = 0.035. The structure of the saddle 
far from the ordered regions is similar in each panel, that of the so-called double fractal sets. The shape of the saddle changes more 
dramatically close to the stability islands. It is evident that the tori do not belong to the saddle. In the numerical simulations A^q = 5- 10^. 
The integration time to used in the algorithm is to = 150 for (a) and to = 180 for {b),(c). 
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Figure 11. Fractal patterns for A = 0.035. (a) Basin boundary. 
A large number A^o = 1-6 ■ 10^ of points was integrated forward 
to see which escape route they chose. The dark gray[in red on- 
line] (white) region represents initial conditions leading to an es- 
cape to minus(plus) infinity. Light gray [orange on-line] marks 
the ordered motion inside the regular islands. The very complex 
structure of the basin boundary indicates chaotic motion. Grid 
size: z ^ 8, Az = 0.008; -2 < 2 ^ 2, Ai = 0.004; (b) Stable 
manifold of the chaotic saddle in the relativistic Sitnikov problem. 
One can see that the filamentary structure of the stable manifold 
and the fractal basin boundary of the escapes are identical as it 
is well-known from dissipative systems. 



set of the chaotic saddle in the phase space. Escape rates 
corresponding to different value of gravitational radii show 
a nearly linear increase with A. The number of trajectories 
with longer life time follows a power-law decay. Such long 
time escape characterizes the so-called sticky orbits which 
may spend very long time in Cantori around the stability is- 
lands. The exponent describing the leakage of sticky trajec- 
tories seems to not depend on A. From a more general point 
of view, our results demonstrate that transient chaos and 
chaotic scattering are robust phenomena in celestial mechan- 



ics as weak relativistic effects are unable to destroy them. 
Perturbations change the characteristic numbers, at most. 
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When deriving the equation of motion, in the correction 
terms we replace the second derivatives with their zeros or- 
der expressions, namely 



kM 
^ = ~2 — r. 
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Further, we use the relation 
1 vr zz 

to obtain 
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APPENDIX A: EQUATION OF MOTION FOR 
THE RSP 

We consider the motion of three gravitating masses in the 
first post-Newtonian approximation of general relativity. 
The first and second particles have the same mass M, op- 
posite velocities, and they revolve in the x — y plane around 
the centre of mass which is at the origin. Their separation 
and relative velocity is denote by r and v, respectively. In 
other terms, 



mi = 7712 = M, ri 



-r2 



^r, vi 



-V2 



(Al) 



The mass 77i of the third particle is negligible compared to 
M. This third particle moves along the z axis. Let us denote 
|ri - ral = |r2 - rsl by p. Clearly, 
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After rearranging terms. 
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In order to obtain the dimensionless equations, first, we 
take the semi-major axis (a) of the Keplerian orbit as unit 
length. Consequently, the unit of the velocity is a/T, where 
T denotes the time unit. As in Section 2.1 we choose T as 
a?'^'^ I (kflMY^"^ corresponding to Kepler's 3rd law (up to a 
constant 27r). Thus, the dimensionless equation from (|A10|) 
is obtained as Equation (|8]) in the main text. 

This paper has been typeset from a T^j]X/ ffl^jX file prepared 
by the author. 



The Lagrangian of the thir d particle can be deriv ed from the 
three-particle Lagrangian fcandau fc Lifsijll975l ) by consid- 
ering r and v as given functions of the time. We get (by 
omitting full time derivatives) 
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The partial derivatives of L are 
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